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A systematic description of multipole degrees of freedom is discussed on the basis of the Stevens' 
operator-equivalent technique. The generalized Stevens' multiplicative factors are derived for all of the 
electric and the magnetic multipoles relevant to /-electron systems. With extensive use of the Stevens' 
factors, we express the spatial dependences of the electric and the magnetic fields, and the electric and the 
magnetic charge densities of localized / electrons. The latter is utilized to draw wave functions including 
their magnetic profile in addition to their shape with the charge density. The definite relation between the 
operators as quantum-mechanical variables in a multipole exchange model and the multipole moments 
in expansion of electromagnetic fields is given. The general treatments for the exchange model with the 
RPA susceptibility and the Ginzburg-Landau free-energy expansion are discussed, using Ce^Lai-^Bj as a 
typical example. The representative formula of the vector spherical harmonics are summarized, which are 
suitable basis for vector fields in the spherical expansion. 

KEYWORDS: multipole, Stevens' multiplicative factor, spherical tensor operator, vector spherical harmonics 



1. Introduction 

Studies with orbital degrees of freedom have encompassed 
a considerable part of condensed matter physics. In fields of 
transition metal oxides a relatively weak spin-orbit coupling 
provides definite description of spin and orbital, and collec- 
tive phenomena are investigated in view of an independent 
or a mutual entanglement of these degrees of freedom. 1 ' 2 ' On 
the other hand, in rare-earth and actinide compounds a strong 
spin-orbit coupling smears spin and orbital forming a harmo- 
nious degrees of freedom called as multipole. 

Because of a localized nature of /-electron wave function 
and its large orbital angular momentum, higher multipoles 
such as octupole become active in orbitally degenerate sys- 
tems. They could play a central role for proper understand- 
ing on mysterious hidden orders with anomalous responses in 
thermodynamics and low-energy excitations. 3 ' 

It is well known that in macroscopic electromagnetism 
the concept of multipole is introduced to characterize source 
charges and currents distributing near the origin. 4 ' An atomic- 
scale counterpart not only characterizes distributions of local- 
ized / electrons but also plays a role of a quantum-mechanical 
variable 5 ' as similar to a spin dipole in ordinary magnetism. 
Such two aspects of multipole degrees of freedom would be 
the cause of possible confusion in a practical calculation. This 
confusion largely arises from indefinite relation between an 
expression in terms of the spherical tensor operator 6 ' and a 
classical definition of the multipole moment. Moreover, there 
exist several different notations with no systematic normaliza- 
tion. A large number of complicated expressions for higher- 
rank multipoles may increase occasional errors as well. 

The purpose of this paper is to give definite and sys- 
tematic description of multipoles and clear relation between 
the quantum-mechanical operators and the classical mul- 
tipole moments in the expansion of electromagnetic field. 
The bridge among these expressions is the Stevens' operator 
equivalent technique 7 ' which is extensively used in an analy- 
sis of energy levels under a crystalline electric field (CEF). 8 ' 
In this paper we generalize the so-called Stevens' multiplica- 



tive factors to all of electric and magnetic multipoles relevant 
to /-electron systems. Using the operator equivalent, we pro- 
vide useful formula for visualization of wave functions in- 
cluding their magnetic profile as well as charge density. In a 
realistic situation with or without uniform external fields, we 
encounter a rather complicated entanglement of plural mul- 
tipoles. To analyze such systems systematically, we discuss 
a general treatment using the random-phase-approximation 
(RPA) susceptibility and the Ginzburg-Landau (GL) free- 
energy expansion for a multipole exchange system. All these 
would be efficient to explore various phenomena concern- 
ing multipole degrees of freedom, and to enhance experimen- 
tal efforts to quantify measurements using NMR, yt/SR, ultra- 
sound, (resonant) X-ray, neutron scattering etc. 

The organization of this paper is as follows. In the next sec- 
tion we demonstrate multipole expansions of scalar and vec- 
tor potentials 9 ' in which we define the classical multipole mo- 
ments. For the given multipole moments we express electric 
and magnetic fields in a simple form with the vector spher- 
ical harmonics. 10 ' The real and the point-group representa- 
tions are also explained. In §3 we discuss the relation be- 
tween the classical multipole moments and the spherical ten- 
sor operators. The operator equivalent technique based on the 
Wigner-Eckart theorem bridges over two expressions yielding 
the generalized Stevens' multiplicative factors. 7 ' 8 ' The useful 
formula for visualization of wave functions are given in §4. 
The fundamental aspect of the multipole exchange system is 
elucidated on the basis of the RPA susceptibility and the GL 
free energy in §5. In §6, we illustrate a use of the present 
arguments with Ce A Lai_ A B6 as an example. The last section 
summarizes the paper. There are three appendices. Appendix 
A contains the definition and representative formula for the 
vector spherical harmonics, which are very useful to express 
vector fields in spherical expansion. The details of the multi- 
pole expansions are given in Appendix B. The derivation of 
the generalized Stevens' multiplicative factors is given in Ap- 
pendix C. To be self-contained and to provide a coherent no- 
tation throughout this paper, we quote several known results 
with appropriate modifications. 
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2. Multipole Expansion 

In this section, we briefly discuss the multipole expansion 
for the scalar and the vector potentials using the spherical 
and the vector spherical harmonics as a basis of the expan- 
sion. Through the expansion, we introduce the electric and the 
magnetic multipole moments, and we examine the symmetry 
property for them. 

2.1 The scalar potential and the electric multipole moment 
Let us start with the Poisson equation for the scalar poten- 
tial, 

V 2 0(r) = -4np(r), (2.1) 

where p(r) is the charge density of localized / electrons. For 
regions outside the source distribution, the solution of the 
Poisson equation is expressed as 



outside the source distribution has the form 



oo p j 

= Z Z — Z P^Qri> W 

P=0 q=-p 

where we have defined the electric multipole moment as 



Qpq = J drr"Z* pq {r)p{r). 
Here, f = r/r is the unit radial vector, and 



Zp q (r) = 



An 
<2p+l 



YUr) 



(23) 



(2.4) 



is the spherical harmonics with the Racah normalization. Note 
that we adopt the Condon-Shortley phase, i.e., [Z pq {f)]* = 
(-\) q Z p - q {r), yielding that Z p o(r) of the odd rank p is a real 
quantity instead of a pure imaginary in the other convention. 
The CEF Hamiltonian with the point-charge model is then ex- 
pressed as 



•HCEF = J] q n <f>(Rn), 



(2.5) 



where q n and R„ represent the charge and the position of the 
ligand ions, respectively. 

The inversion operation transforms p{r) to p(-r). Using 
Zpq(-r) = {-lyZpqif), we show that Q pq is transformed to 
{-^YQpq- Namely, the electric multipole moment has the par- 
ity (-l) p . If a system has the inversion symmetry, the odd- 
rank electric multipole moments vanish since p(r) = p(-r) 
holds. The time-reversal operation changes nothing on p(r). 
Therefore, the electric multipole moment is even under time 
reversal. Since p(r) is real, we have Q* pq = {-l) q Q p - q . 

2.2 The vector potential and the magnetic multipole moment 
The Poisson equation for the vector potential is given by 



V 2 A(r) = j(r), 



(2.6) 



where j(r) is the current density originating from the orbital 
and the spin currents of / electrons. In contrast to the scalar 
potential, the vector potential has an intrinsic angular momen- 
tum ("spin") 1 . Taking this property under consideration, we 
require basis vector fields which transform like Z pq (r). One 
suitable basis is known as the vector spherical harmonics. 10) 
In the gauge V • A = 0, the expansion of the vector potential 



^ A 1 [ tZ pq (r) \ , 



(2.7) 



P=0 q=-p 

where the magnetic multipole moment is given by 

M pq = j drV[r<>z; q (r)] ■ M(r). (2.8) 

Here, t — —ir x V is the (dimensionless) orbital angular mo- 
mentum, and (Z pq (r) is one of the vector spherical harmonics 
(without normalization) as shown in Appendix A. M(r) de- 
notes the magnetization density defined through 



j(r) = cV x M(r). 
By the partial integration in (2.8), we have 



(2.9) 



pq - j ' drr"z; q (r)p m {r), (2.10) 

where we have introduced the magnetic charge density, 10) 

Pm (r) = -V • M(r). (2.11) 

This expression is formally similar to that of the electric mul- 
tipole moment, (2.3). 

In contrast to the charge density, the magnetic charge den- 
sity as well as M pq is odd under time reversal. The inver- 
sion operation transforms M{r) to M(-r), and consequently 
Pm(r) -> -p m (-r). Thus, the magnetic multipole moment has 
the parity (-l) p+l . In the presence of the inversion symmetry, 
the even-rank magnetic multipole moments also vanish since 
M(r) = M(-r) [p m (r) = -p m (-r)]. The complex conjugation 
isM; e = (-l)«M p _ ? . 

2.3 The electric and magnetic fields 

From (A- 13d), we have the simple expression of the electric 
field in the multipole expansion, 



p=0 q=-p 

(2.12) 

Similarly, we obtain the expression of the magnetic field, 



co p 



JKr) = V x Air) = -2 J 2 J 

p=0 q=-p 

(2.13) 

which is formally similar to E(r). For the given electric and 
magnetic multipole moments, the electric and the magnetic 
fields are calculated straightforwardly by using the definition 
of the vector spherical harmonics, (A-l). It is also useful to 
give the scalar and the vector products with r, 

>-E{r) = f i f i (p + l)Q„^, (2.14) 

p=0 q=-p 

'**w = ZI> 9 ^. (2.15) 

p=Q q=-p 



Similar expressions for B(r) are obtained by replacing Q 



pq 



with M, 



pq- 
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2.4 The real and the point-group representations 

The multipole with q + in the spherical representation is 
a complex quantity. It is useful to introduce the real represen- 
tation for q > as follows: 



4(c) 
A pq 



V2 



AM - (a _a*\ 



(2.16) 



where A pq represents any quantity that transforms like Z pq {r) 
under spatial rotation. The corresponding real expressions for 
r p Z pq (f) (called as the tesseral harmonics 8 ') are summarized 
in Table C- 1 . In the real representation, a sum of products is 
rewritten as 



X A pi B p 



q=-p 



£ [Af q Bf q + Af q Bf q ] + ApoBpQ. (2.17) 

q=\ 



A similar transformation is applicable to vector fields as well. 

In reality, magnetic ions are placed in a crystal with a 
proper point-group symmetry. When a CEF splitting is small 
and a total-angular momentum (/) multiplet can be treated 
as a whole, the spherical or the real representation is ap- 
propriate. On the other hand, when a CEF splitting is large 
and one of CEF multiplets dominates low-energy physics, the 
point-group irreducible representation is suitable to classify 
the multipole moments. Due to the fact that any point group is 
a subgroup of the rotation group, the point-group harmonics 
are constructed as linear combinations of the spherical har- 
monics. For instance, the cubic harmonics 11 ' 12 ' under Oh are 
given in Table C-2, where F and y represent the irreducible 
representation and its component in the Bethe notation. 6 ' Note 
that any scalar product in the spherical expansion can be re- 
placed by the point-group identity representation, which has 
the form with a sum of pairs of the same irreducible represen- 
tation, 



r 

2 

q=-p 



A* R 



Ty 



plyBpTy 



(2.18) 



In the following sections, we often give results only in the 
spherical representation. Any scalar product may be replaced 
properly by the corresponding point-group representation. 

3. Multipole and Stevens' Operators 

3.1 Multipole operators 

In the previous section, we have introduced the electric and 
the magnetic multipole moments, which are determined by 
the charge density p(r) and the current density j{r) [or equiva- 
lently the magnetization density M(r)] of the localized / elec- 
trons. In the quantum statistical mechanics, the charge and the 
current densities should be regarded as a thermal average over 
/-electron states. 

The corresponding one-body charge density operator act- 
ing on /-electron wave functions is given by 



p(r) = -ej]6(r-r j ), (e > 0), 



(3.1) 



where the summation is taken over all / electrons. With this 
operator, it is natural to introduce the electric multipole oper- 



ator as 



Qp, = ~e X J dr 6(r - rj^Z^rj). (3. 



2) 



A derivation of the magnetic multipole operator is more 
involved. The detailed discussion is left in Appendix B, and 
we quote the result, 

2(j 



M„ 



p + l 



2s, 



(3.3) 

where p B = -eh/lmc is the Bohr magneton, and (j and Sj are 
the orbital and the spin operators of j-th / electron. 

With these operators the classical multipole moments are 
given by the thermal average over /-electron states, 



Qpq - (Qpq) , > 



M 



= (M pq ) f . 



(3.4) 



3.2 Spherical tensor and Stevens ' operators 

In order to calculate systematically a matrix element of the 
multipole operators, let us consider the spherical (Racah) ten- 
sor operator, 6 ' J pq , which is defined by the /?-th polynomial 
of the total angular momentum operator, / = (J x , J y , J z ). The 
definition of the spherical tensor operator is 



f PP = (-!)" ■ 



(2p-l)!! 



(2p) ! 



(j + y, 



[/_, J pq \ = yj{p + q)(p -q+ l)J pq -u (q<p), (3.5) 

where J ± = J x ± iJ y . We express the Wigner-Eckart theorem 
for the spherical tensor operator, 



(jM'\j pq \jM) = {-\) J+M -''\ [ _ J MI 



J 
M 



gwii'>. 

(3.6) 

where the parenthesis denotes the 3j symbol, 13 ' and the re- 
duced matrix element of J pq is given by 



\ II p\\ I 2 p\ (2J-p)\ 



(3.7) 



With use of the Wigner-Eckart theorem, we compute any ma- 
trix element of J pq within a J multiplet. 

It is also possible to construct the hermite tensor operator 8 ' 
from J pq with a similar linear combination as (2.16). The ex- 
plicit expression of the operator can be obtained by replacing 
(x, y, z) in r p Z ( p c q , r p Z p ^ q and r p Z p o with the symmetrized prod- 
uct of (J x , J y , J z ), i.e., 

x k y m z n 



k\m\n\ 



(3.8) 



(k + m + ri)\ ^ 

where the summation is taken over possible permutations. For 
instance, we obtain the operator form of by replacing 



r 3 Z (c) = 



VB 2 2 
-^-20 - y ) 

Vis 



Jy) (J x Jy)Jz J X-J X JyJ "<A j • 

(3.9) 



The point-group counterparts are obtained in a similar way. 
We compute easily any matrix elements of J pq and its variant 
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with this prescription. 

The so-called Stevens' operators 7 ' 8 ' are often used in the 
CEF theory. They are proportional to the hermite tensor oper- 
ators as follows: 



o 2 ° = 2j 20 , d 2 2 



2 «c) 

V3 22 ' 



0\ - 8/40, - —^J^2' 



4 



J 43 ' 



04 (s) = 



m 

- -/ 43' 



V70 

06° = 16^60, Ob 2 = 
16 



f(c) 

V35 44 ' 



V70 



32 



. tic) 

3V7 64 ' 



/210 
32 



f(c) 
J 62' 



16 



/210 



V35 

f(c) 
'63' 



f(s) 
''44' 



/462 



(3.10) 



Note that the spherical and the point-group tensor operators 
are properly normalized, but the Stevens' operators are not. 
The former is more appropriate for systematic calculation. 

3.3 Generalization of Stevens ' operator equivalents 

Since J ' pq is the p-th polynomial of the axial vector J a with 
time-reversal odd, the spherical tensor operator, J pq , has the 
even parity and the time reversal (-1)'' Hereafter, the time- 
reversal symmetry of the tensor operator will be indicated 
by g(erade) and u(ngerade) in the subscript. Moreover, the 
spherical tensor operator is transformed as Z pq (f) by defini- 
tion. These symmetry properties are common with the even- 
rank electric and the odd-rank magnetic multipole operators. 
Note that the even parity of the tensor operator differs from 
(— l) p parity of the spherical harmonics. Therefore, according 
to the Wigner-Eckart theorem, we conclude that any matrix 
elements of the multipole operators are proportional to those 
of the corresponding spherical tensor operators. We assume 
the presence of the inversion symmetry in what follows, then 
all of the relevant multipole operators are described by the 
spherical tensor operators. 

Let us express any matrix element within a J multiplet in 
an /" configuration as 

^{jM'\j pq \jM), (3.11a) 



(nJM' | Q pq | nJM) 



(nJM 1 1 M pq | nJM) = p B (r^ 1 ) g { n p) (jM' \ J p 



\jm). 



where we have defined the radial average, 



(3.11b) 



(3.12) 



The relativistic Hartree-Fock estimate for (r k ^ may be found 

in ref. [14], for example. Here, g„ is the generalized Stevens' 
multiplicative factor, which is independent of M, M' and q. 
For the Hund's-rule ground multiplet, J L 2 s+i, in the Russell- 
Sanders (LS ) scheme, g„ is given in Table C-3, and shown in 
Figs.l and 2 as a function of the number of / electrons. The 
ratio of the orbital and the spin contributions to the magnetic 
multipole, 



(p) _ gn ] (orbital) 

n ~ (P) ' 

8n (spin) 



(3.13) 
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Fig. 1. (Color online) The even-rank Stevens' factors as a function of n. 
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Fig. 2. (Color online) The odd-rank Stevens' factors as a function of n. 



is shown in Fig. 3. The derivation of the generalized Stevens' 
factor is given in Appendix C. 

Note that the even-rank g^ is equivalent to the ordinary 
Stevens' factor, 8 ' and gj, 1 ' is nothing but the Lande's g fac- 
tor, 6 ' 8 ' i.e., 



= gj, 



(3.14) 



and 



£1 = e 2 =aj, &> = 0i=fij, g^ = e 6 =yj. (3.15) 

Since the scalar and the vector potentials are one-body 
fields, the multipole operators interacting with them are rep- 
resented by one-body operators as shown in (3.2) and (3.3). 
Due to the selection rule (CT7) of the spherical harmonics, 
the Stevens' factors more than rank p > 21 - 6 must vanish. 
In view of a quantum-mechanical transition, multipole oper- 
ators with p > 6, e.g., <9 8 8 ~ (|4, +4) (4, -4] + h.c), could 
exist for J > 3. However, such operators with p > 6 are rep- 
resented by more-than two-body operators, and do not couple 
with one-body potentials. 

Let us express an arbitrary state within a J multiplet, 



\y) = Yj u vy l /M > ■ 



(3.16) 



Note that the unitary matrix Umj could be complex. With the 
operator-equivalent method, the classical multipole moments 
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Fig. 3. (Color online) The ratio of the orbital and the spin contributions to 
the magnetic multipoles. 



M " = ^(/||^)£^ ? (y)-^ 

y 



are given by 

Qm 

-e(rP)' ^(rP-i) °" m-ir/^ Zf 

(3.17) 

where we have defined the weight function for a y state as 
W„(y) = £(-l)' + *-'(_£, M P } U M'yUMy (3.18) 

MM' ' ^' 

Note that in the case of U My - S My , *Wpq(y) = unless q = 0. 

4. Visualization of Wave Function 

In this section, we consider visualization for the charge 
density p(r;y) and the magnetic charge density p m (r;y) of 
/ electrons in a particular state \y). From (2.3) and (3.4), we 
may obtain the relation 

(y | Q Pq \j) = fdr rPZ; q {r)p{r- y). (4. 1 ) 

Suppose that we write p(r; y) in the separable form, 

p(r; y) = -eR 2 f (r)p e (r; y)/4n. (4.2) 

Substituting this into (4.1) and using the completeness of the 
spherical harmonics, 

2p + 1 , 



z 

PI 



z pq(r)Z*(r) = 5(r-f'), 



4 n ^PlV >^pq 



we obtain the angle dependence of the charge density, 

DO P 



s--, x , (y | Q Pq I y) 

p e (r;y) = J] J] (2p + l) -'-^ Zp « {r) - 

p=0 q=-p 



-e (rP) 



(4.3) 



(4.4) 



This expression is equivalent to that obtained by Walter. 15 ^ 
Similarly, supposing that 



R)(r) 

Pm(r; y) = pb Pm(h y)/4n 



(4.5) 



and (2.10), we obtain the angle dependence of the magnetic 
charge density, 



oo p 



Pm(r;r) = 2Z (2 ^ +1) 

P=0 q=~P 



(y\M pq \y) 



ZnJr). (4.6) 



M^- 1 ) pqy 

Using the operator equivalents (3.1 1), we express the angle 



dependences in the common form, 

6 p 

Pe, m (r; r) = J](2 P + Dgf (J || J P || J) 2 nV M (y)Z M (r), 

p=0 q=-p 

(4.7) 

where the summation for p is taken over even (odd) integers 
for p e (p m ). This expression is useful to visualize a wave func- 
tion for a y state. The presentation of the charge density is 
similar to that in refs. [15, 16]. Namely, the radius R y (r) to the 
surface of the 3-dimensional plot is defined as 



R y {f) = [p e (r\y)f , 



(4.8) 



in which a — 1 is chosen to emphasize the gradation of the 
charge density, although a = 1/3 is natural to yield that the 
encircled volume becomes the total charge. The surface-color 
map is used to represent the magnetic charge density p m (r; y), 
which is normalized to hold the range of the distribution in 
[-1,1]. 

5. Treatments for Multipole Exchange Systems 

5.1 The exchange model in a crystal 

As was mentioned in §2.4, the point-group representation 
is appropriate in the case of large CEF splitting. When one 
of CEF multiplets or a bunch of CEF multiplets with small 
splittings dominate low-energy physics, we consider a multi- 
pole exchange system within the relevant CEF states. In the 
restricted basis with the (pseudo) degeneracy d, the multipole 
operators in the point-group representation become reducible, 
and some of them are proportional with each other. 

Meanwhile, in the view of the quantum-mechanical vari- 
ables, we require d 2 independent operators to expand the re- 
stricted manifold (one of these is the identity operator). Since 
d 2 = 2^=o(2p + 1), we formally assign the spherical tensor 
operators up to rank d - 1 to the d 2 independent operators. 5) 
Thus, such mathematically independent operator is called the 
multipole operator as well. 

Keeping this consideration in mind, we denote the d 2 - 1 
independent operators at the site i as X? {a = 1, 2, • • • d 2 - 1), 
except the identity operator denoted as 1, . The operator X" is 
hermite, traceless and is normalized as 



(5.1) 



In terms of X a , we write down a generalized exchange model 
with uniform external fields, 



Mr 



ij a/3 ia 



where we assume complete degeneracy of CEF states for sim- 
plicity, but a generalization to pseudo-degeneracy is straight- 
forward. Note that possible differences in normalization of the 
multipole operators are absorbed in the definition of the ex- 
change coupling Df: . However, it should be emphasized that 
the relation between X" and the multipole operators, Q pq and 
M pq , is vital in evaluating an effect of the multipoles through 
electromagnetic probes. The non-vanishing combination of a 
and p in may be obtained by symmetry consideration with 
respect to an interacting bond i-j. l7) The second term repre- 
sents coupling with uniform external fields such as magnetic 
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field, uniaxial strain, and so on. 

As was mentioned in §3.3, the multipole operators with 
p > 6 consist of more than two-body operators. The ex- 
change couplings are expected to be small for such opera- 
tors because the origin of Df? involves higher-order exchange 
processes to transfer more than two electron states simulta- 
neously. To the contrary, Dff for p < 6 could be the same 
order in magnitude. This is because the origin of the coupling 
is the RKKY and/or the superexchange mechanism in which 
the matrix elements and the intermediate energies have simi- 
lar strength. Especially, when a virtual process through a fea- 
tureless state such as f° configuration dominates, all of the 
coupling strength turns to be coincident with each other. 18 ' 

Since the operators Xf span the restricted manifold, the 
product of the operators can be expanded as 

W = Z + + 8a ^ (53) 

7 

where the symmetric and anti-symmetric structure constants 
are calculated from the definition of Xf as 

8a* = ^Tr ; ([xfxf + Xf Xf] XJ) , (5.4a) 

ifapy = ^Tr, ([Xf*f - xfxf ] XJ) . (5.4b) 

The symmetry property of the system is completely deter- 
mined by the structure constants. In the case of the Pauli ma- 
trices, X" = &f, we have g a p 7 = and f aj s y = e aj3y , where e a p y 
is the anti-symmetric (Levi-Civita) symbol. 

5.2 The RPA susceptibility 

Let us consider the static susceptibility of the multipoles 
within RPA. The second-order phase transition from a disor- 
der phase is then determined by the divergence of the suscep- 
tibility in the mean-field approximation. For this purpose, we 
add a coupling with fictitious fields to the exchange Hamilto- 
nian, 



(5.5) 



We divide the thermal average of the multipole operator into 
two parts, 



(5.6) 



where the second term is the induced moment being propor- 
tional to the infinite simally small field <f> a . Note that (Xf ^ 
could be finite in the presence of external fields, otherwise 
they must vanish in the disorder phase. 
The mean-field Hamiltonian is 



•Wmf = X <H f (i)-J]xfAf 



(5.7) 



where we define the mean-field Hamiltonian without <pf 
fields, 



and the effective infinitesimal field 

P 



Xf, (5.8) 



(5.9) 



The self-consistent equation for h a + and cpf = is given 
by 



{ jt ?) f = Y J M m \ jt ?\ m )> 



2m ^ 



-0E„ ' 



(5.10) 



where 'Hf(i) \m) = E m \m) and \Xfj^ is independent of the 
site in the disorder phase. 

According to the linear-response theory, we have the local 
susceptibility for <p" = 0, 

^=Z/- 1 V^" J " ("i^i")("l^l") 



Note that in the absence of the external fields, h a = 0, we have 
the Curie law,^ c (/z ff = 0) = fi6 aj3 . Using we obtain 



B y V kS 



(5.12) 



By the definition of the susceptibility, yf. - d (Xf) /<9<^L=o, 

l J \ 1 I (j) J 

we have the relation, 



kS 



(5.13) 



With the Fourier transformation, we finally obtain the RPA 
susceptibility, 



yP 



(5.14) 



where D a P(q) = e _,?r "Z) a ?. In the case of diagonal cou- 
pling D a P = D a 8 a p and h a = 0, we have the simple Curie- 
Weiss susceptibility, 



T - D a (q) 

The second-order phase transition is determined by 



det 



<W-IXc^(<7) 



0, 



(5.15) 



(5.16) 



at q — Q with the maximum T c . The ratio of the order param- 
eters just below T c is determined by the eigenvector of the 
matrix in the eigenvalue equation. 

5.3 The Ginzburg-Landaufree energy 

The GL free-energy expansion is useful to elucidate sys- 
tematically an entanglement of the multipoles. Here we derive 
the general expression of the GL free energy of the exchange 
model for h a = 0. Let us consider the one-body trial Hamilto- 
nian, 



(5.17) 



Then, the exact free energy is upper-bounded by the Feynman 
inequality, 19 ' 



T<T tI = To + {'Hf-'Ho) , 



(5.18) 



where (• • -) is the thermal average with respect to 'Ho and 
To = -/T 1 In Tre^ . Minimizing Tu with respect to ifff, we 
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obtain the best Hamiltonian within the one-body approxima- 
tion. This type of variational treatment is equivalent to the 
mean-field theory. 

The thermal average is easily evaluated as 

ij a/3 ia 

(5.19) 

By using d(X"J Q /dt//l = f36ikd ay , the stationary condition 
yields 



&Fbh/dX" = 0, and it is given in the GL expansion as 



p 



+ ^J]L aP7 XfxJXf = 0. (5.27) 



PyS 



*? = Z«*»c 



The fluctuation from the stationary is related to the suscepti- 
bility. Namely, we replace X" — > X" + 8X" and retain terms 
up to 2nd order in 6X", then we obtain the deviation of the 
free energy, 



<A?=<A, 



(5.20) 



which plays a role of the self-consistent equation. Eliminat- 



(5.28) 



ing the order parameter X" = (xf) Q _^ with the stationary (*" % ~ TS ijd<*p 



condition, we obtain the best trial free energy as 

MzzmT^Z^K 2 ^ 

n,ft " i 



(5.21) 



Now, we express it in terms of the order parameter X". Ex- 

—3 

panding X" up to \p with use of (5.3), we have 



-rP-rr 



+ %■ Z ( L af3y5 ~ ISa/jSyd) iffi^, (5.22) 



I3y5 

where we have introduced 



L a f3y5 - ^ gaftgydZ + S a p8- 



yd' 



(5.23) 



The converse relation of (5.22) is obtained by a recursive ex- 
pansion as 

0] ~ Xf - \ £ ga/iy XfXj + 1 2 L^^Xf. (5.24) 



Using this relation and (5.20), we finally obtain the GL free 
energy up to 4th order (except the -TN In d term), 



Xf 



1/038 



iafiyd 



In the case of X? = cr" and = Jij5 a fj, we recover the 
GL free energy of the SU(2) Heisenberg model (we denote 

< = x f)> 

1 7" 
^gl = ^ Z ~ ^) ' m J + 12 Z (nl < ' ,Mi)2 - (5 ' 26) 

The self-consistent equation (5.20) reduces to the condition 



Z 8a/3yXj - T Z ( Lff ^ + 2L «y/») ^T^f 



L 7 



(5.29) 



Since X" = in the disorder phase, we recover the RPA sus- 
ceptibility (5.14) for h a = 0. On the other hand, the mode 
mixing arises through the 3rd and 4th-order couplings in the 
ordered phase, X" + 0. 

The simplest non-trivial entanglement of the multipoles 
arise from the 3rd-order coupling, g a p y . The coupling con- 
serves the momenta of the order parameters. Namely, the en- 
tanglement comes from uniqueness of the local symmetry- 
breaking wave function. When a spontaneous order occurs for 
Xf, the mode mixing takes place between xf and Xj. If two 
of three multipoles are equivalent, e.g. a = f3, Xj is induced 
eventually. Since the free energy is even under time reversal, 
the 3rd-order term should consist of three electric multipoles 
or two magnetic and one electric multipoles. Consequently, 
the magnetic multipole is the primary order parameter while 
the electric multipole is secondary when both multipoles co- 
exist in the ordered phase. 

6. The example with Ce^Lai-^Bfi 

In this section, we illustrate a use of the previous arguments 
with Cei_ A La t B6 as an example. The trivalent Ce ion is placed 
in the cubic crystal field Oh, and J = 5/2 multiplet splits into 
T 7 doublet and Tg quartet. The wave functions are given by 



= V¥ hi) - Vf 



ir 7 ;±> 



|r 8 ;a+> = -J- 



4WS 



\n-,b±) 



1 

+ — 

2 



(6.1a) 



(6.1b) 



(6.1c) 



T 7 doublet consists of the time-reversal pair, while r& quartet 
has an additional orbital degrees of freedom. The splitting be- 
tween Fg ground state and F 7 excited state is about 500K, so 
that it provides an ideal quartet system in low temperatures. 
There are four phases at most in H-T phase diagram, which 
are called as I-IV. 20) 

Within T 7 doublet, we decompose the direct product of the 
basis as F 7 ® T 7 = Fi e T 4 , which indicates that two types of 
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the multipole operators including the identity one are active 
in F 7 doublet. We express a set of independent operators in 
terms of the 2x2 Pauli matrices <x y and the identity matrix 1. 
They are related to the spherical tensor operators as 

T\g : 1 = Jo,ig,i = —J. 



5 V21 



4,lf,l> 



(6.2a) 



X 



1,2,3 



1 



— °"y — ^ J 1 ,4u,y — Jq J 3,4u,y 



65 15 V35 



J5,4hu,y, (6.2b) 



where the cubic tensor operator is expressed as f P ,r(g,u),y Note 
that the operators belonging to the same F with different ranks 
are reducible in the subspace. 

On the other hand, the decomposition for F 8 quartet 5) is 
r 8 ® T 8 = Fi e r 2 © F 3 © 2F 4 © 2F 5 . A set of independent 
operators in F 8 quartet is given by 

- t 2 



2u 



5 V2T 



X 1 = ty = 



F 3g : X 2 ' 3 = (t z ,t x ) = T^2,3g,y = - 



9 V5 



■^3,2«,1. 



1 



r . V 4 > 5 > 6 — fir 



4 

14 .. 

~5 J 



5 Vl5 



■^4,3g,: 



(6.3a) 
(6.3b) 
(6.3c) 



l,4u,y ' 



45 



4u,y 



— J5,4au,y ' 



T4BU '■ x 



7,8,9 



(fj+& x , f)-& y , t z fr z ) 



9 V35 
2 



'5,4bu,y> 



(6.3d) 



15 



J\,4u,y + TT-^3 



45 



4u,y 



2 f 14 , 

75 45 V35 

F 5 „ : X mi ' 12 = {l + & x ,l-& y ,f x & z ) 



(6.3e) 



Table I. The symmetric structure constants, g a py 





Saf3y 


a,P,y 


Safiy 


a,P,y 


Safiy 


1,4,13 


1 


1,5,14 


1 


1,6,15 


1 


2,4,7 


-1/2 


2,4,10 


V3/2 


2,5,8 


1/2 


2,5,11 


-V3/2 


2,6,9 


1 


3,4,7 


-V3/2 


3,4,10 


-1/2 


3,5,8 


- V3/2 


3,5,11 


-1/2 


3,6,12 


1 


7,8,15 


-V3/2 


7,9,14 


V3/2 


7,11,15 


1/2 


7,12,14 


-1/2 


8,9,13 


-V3/2 


8,10,15 


1/2 


8,12,13 


-1/2 


9,10,14 


1/2 


9,11,13 


"1/2 


10,11,15 


V3/2 


10,12,14 


V3/2 


11,12,13 


V3/2 










Table 


II. The anti-symmetric structure constants, f c 


r/3y- 


a 8 y 


fapy 


a,p,y 


fafiy 


a,fi,y 


fctfty 


1,2,3 


1 


1,7,10 


1 


1,8,11 


-1 


1,9,12 


1 


2,7,13 


-V3/2 


2,8,14 


-V3/2 


2,10,13 


-1/2 


2,11,14 


-1/2 


2,12,15 


1 


3,7,13 


1/2 


3,8,14 


-1/2 


3,9,15 


-1 


3,10,13 


-V3/2 


3,11,14 


V3/2 


4,5,6 


1 


4,8,9 


1/2 


4,8,12 


-V3/2 


4,9,11 


V3/2 


4,11,12 


-1/2 


4,14,15 


1 


5,7,9 


1/2 


5,7,12 


V3/2 


5,9,10 


V3/2 


5,10,12 


1/2 


5,13,15 


-1 


6,7,8 


1/2 


6,7,11 


V3/2 


6,8,10 


-V3/2 


6,10,11 


-1/2 


6,13,14 


1 




r 



5g 



X 



13,14,15 



3 V5 



J%5u,y 



15 V35 
1 

5 Vl5 



Js,5u,y, 
■^4,5j,y> 



(6.3f) 



(6.3g) 



where & y and f 7 are 4 x4 matrices acting on the time-reversal 
and the orbital pairs respectively, and we have defined 



fj ± = ~(t x ± V3r z ), 



(6.4a) 



(6.4b) 



Note that the magnetic dipole moment belongs to irre- 
ducible representation, and is expressed as 

4 



My 



„(1) J _ \'4..XO . ■ y 

8i J\M,y ~ A 7 



7,8,9 



(6.5) 



The operators X a satisfy the orthonormality, 2I) (5.1). The 
symmetric and the anti-symmetric structure constants are 
summarized in Tables I and II. 

From (4.7), the angle dependences of r 7 wave functions are 
expressed in terms of the cubic harmonics, Z p f, 7 (r), 



p e (r;F 7 ±) = gf> - 45 ^2Agfz AXl {f\ 



Fig. 4. (Color online) The charge and the magnetic charge densities for 
Ce 3+ viewed from [111], (a) T 7 ±, (b) r 8fl ± and (c) r 8i ,±. 



p m (r;F 7 ±) = + 



^ ( 1 1) Z 1 , 4 , 3 (r) - 70^ 3, Z3, 4 ,3(r) 



~gf } {l3Z 5 , 4fl , 3 (r) - 3 V35Z 5 , 4W (r)} 



(6.6b) 



Since r 7 doublet has no orbital degrees of freedom, p e (r; r 7 ±) 
is expressed only with the even-rank Fj harmonics. On the 
other hand, p m (r; F 7 ±) consists of the odd-rank F 4 harmonics, 
which belong to the same irreducible representation of the 
magnetic dipole moment. Taking an average over the time- 
reversal pair, we have no magnetic charge density for F 7 
states. F 7 wave functions are shown in Fig. 4(a). 

Similarly, F 8 wave functions are shown in Fig. 4(b) and (c). 
Taking an average over F 8 quartet, we have the full-symmetric 
charge density, 

Pe(r;r 8 ) = l - J]p e (r;r 8r ) = gf + y V2T^ 4) Z 4 , M (r), 

y 

(6.7) 

(6.6a) an( j jjj e magnetic charge density vanishes. The averaged wave 
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Fig. 5. (Color online) The charge and the magnetic charge densities of 
in the ordered phase, (a) disorder, (b) Tj g ^, (c) T2,, and (d) rj." 1 '. 



function is shown in Fig. 5(a). 

The high-temperature phase II is considered as the antiferro 
(AF) r*, g quadrupole order. 22 ' Taking z axis as a quantization 
axis, r5 g 3 molecular field lifts the quartet into two doublets. 
The lower doublet is expressed as 

\5g3; ±> = -L [± |F 8 ; a+) - i |F 8 ; b±}] . (6.8) 
V2 

An average over the doublet states yields the charge density, 



1 



p e (r; 5g3) = | ^ p e (r; 5g3y) = p e (r; F 8 ) 

y 

- 5gfz 2 ^(r) - 45 VT5g ( 1 4) Z 4 ,5,3(r), (6.9) 

while the magnetic charge density vanishes. It is natural that 
the deviation from the cubic symmetry is characterized by the 
even-rank Fs y — 3 harmonics. The averaged wave function 
is shown in Fig. 5(b). 

In Table I, there is the third-order coupling among 2m- 
4Au3-5g3 multipoles in the GL free energy. When we apply 
a magnetic field in [001], the uniform magnetic dipole mo- 
ment (T^ijfl) arises. As a result, the AF F2„ octupole is induced 
through the 3rd-order coupling. Thus, the phase II is more sta- 
bilized as the magnetic field increases, provided that the AF 
2u-2u exchange coupling presents. Similarly, a magnetic field 
in [110] with the AFr5 g i+r5 g 2 order induces the AF magnetic 
dipole moment (T4A«3), which was observed by the neutron 
scattering. This mechanism was discussed by Shiina et al. in 
the mean-field approximation, 5 ' which first indicated that the 
magnetic octupole plays an important role behind anomalous 
phenomena. 23 ' The low-temperature phase III 24 ' was exam- 
ined similarly by the extensive use of the GL expansion ac- 
companied with r*, g AF quadrupole ordering. 25 ' 26 ' 

We briefly mention the pure Ti u octupole order, 27 ' which 
has not been observed so far. F2« molecular field also lifts the 
quartet into two doublets. The lower doublet is given by 



\2u;±) = -^-Un;a±)-i\T s ;b±)] 
V2 



(6.10) 



This state breaks the time-reversal symmetry, however, the 
uniform magnetic susceptibility remains increasing with de- 
crease of temperature in contrast to the ordinary cusp-like be- 
havior of a magnetic order, since there still exist two-fold de- 
generacy. We have the averaged charge and magnetic charge 
densities, 

p e (r;2u) = )-^p e (r;2uy) = p e (r;F 8 ) (6.11a) 



p m (r;2w) = - ^p m (r;2uy) = -— V5gf ) Z 3i 2,i(?'). 
y 

(6.11b) 

The charge density has the full crystal symmetry, while the 
magnetic charge density is characterized by the odd-rank F2 
harmonics. The wave function is shown in Fig. 5(c). 

Finally, we consider the phase IV, which is considered as 
the AF F^,, magnetic octupole phase. 25 ' 2S ~ 32 ' When Fs u molec- 
ular field is applied along the high-symmetry axes, [001], 
[110] and [111], the maximum eigenvalue is obtained in 
[111]. Namely, the easy axis for Fs u magnetic octupole is 
[111]. Thus, we consider AF Fs,, magnetic octupole order, in 
which the operator, 

1 : (x, 10 +X, 11 +X, 12 ) (6.12) 



X?" = 



V3 



becomes diagonal with singlet-double-singlet eigenvalues. 
The non-degenerate ground state gives rise to the cusp-like 
behavior in the uniform magnetic susceptibility. 20 ' The low- 
est singlet is given by 

|r 8 ;a+> + — |F 8 ;a-> 



|5m) = 



1 



(7+4 V2!') 1 



V3 



V3 



-(-l) 1/4 |F 8 ;Z?+> + |F 8 ; /?->], (6.13) 
and both densities are expressed as 

p e (r; 5m) = p e (r; F 8 ) - 5gf%^r) - 45 Vl5g^Z 4 , 5 (r), 

(6.14a) 



p m (r; 5m) = -21 VTo/ 1 3) Z 3 , 5 (r) - ^ V70^f ) Z 5; 5(r), 



(6.14b) 



where 



Z pX (r) = — [Z p , r ,i (r) + Z pXa (r) + Z pX ^r)\ . (6.15) 

The deviations from the cubic symmetry are characterized by 
the even-rank Fs harmonics and the odd-rank harmonics, 
respectively. The wave function is shown in Fig. 5(d). 

Let us denote the primary order parameter, <p = X^, U {Q) = 
[X 5u i(Q) + X 5u2 (Q) + X 5u3 (Q)y V3. There is the coupling be- 
tween Tsj, magnetic octupole and the r$ g electric quadrupole 
with the principal axis [111]. We denote the secondary or- 
der parameter as £ = X 5g (0). We only consider the nearest 
neighbor AF couplings D"".^ = -J a < for 5m-5m and 5g-5g 
multipoles. Then, the relevant GL free energy is given by 



a 



Tcl = ~(T - 7»0 2 + ^ 4 + a -f + c<p 2 % + ■■■ 



(6.16) 



where we have defined the critical temperatures, T a = 6J a . 
We introduce parameters, a, a, b and c, which are evaluated 
from the symmetric structure constants as a — 1, b ~ 27^/3, 
a ~ T^ + Tg, c ~ T<f,/2. Minimizing the free energy, we have T 
dependences of the primary and the secondary order parame- 
ters as 

<m = ^A(T+ - T), f(D = -°- [<P(T)] 2 , (6.17) 

where A = aa/(ab - 2c 1 ) > 0. The primary order parame- 
ter has the ordinary square-root T dependence, while the sec- 
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ondary has T-linear dependence. 29 ' The induced uniform or- 
der £ gives rise to the lattice distortion in [111], which is in- 
deed observed experimentally. 33 ' 

The contribution from the fluctuations to the free energy is 



with the inverse matrix of the susceptibility, 

(a(T -r ) + 3/?0 2 + 2c£ 2c(f> 
2c(f> a 



r 1 



(6.18) 



(6.19) 



Then, the susceptibilities in the disorder phase are given by 
X*(T)= * - - , Xf(T) = -, (J>T 4> ). (6.20) 



In the ordered phase, we obtain the susceptibilities as 

1 bA 

, X( (T) = — , 

laiTs - T) aa 



(r<7». (6.21) 



The susceptibility of the primary order parameter is divergent 
toward T$. On the other hand, of the secondary order pa- 
rameter has a discontinuity at Tq, since the correlation length 
for £ remains finite. The discontinuity is given by 

2c 2 



4% =*f(Z>-)-Ar f (7V+) 



rA > 0. 



aa'- 



(6.22) 



Note that vanishes when c — > 0. The change of the elastic 
constant in C44 mode corresponds to ~x%. Therefore, the pos- 
itive jump ofxt leads to a sudden softening of C44 mode. 34 ' 

The magnetic octupole order yields the internal magnetic 
fields around the Ce ion. 35 ' 36 ' From (2.13), we obtain the in- 
ternal magnetic field, 37 ' 



■M 5 , 5u F 5 6 5 „(r), (6.23) 



B{r) = --^M X5u Yi 5u (r)- 

where the magnetic multipole moments are expressed from 
(3.11)and(6.3f)as 



M 3 , 5 „ =yUB (r 2 )^ 3) -3V50, 



(6.24a) 



M 5>5 „ = fi B {r A )g { [ 



(5) 



15 V35 



0, 



(6.24b) 



and the relevant vector spherical harmonics are given by 



Yi 5u (r) 



z « " Vf Z4 °) cx " z « " Vf Z4 °) 



'41 



(6.25a) 



V5 



V2T "1 f (s) V21 , 

— Z 60 e, - + ~^-Z 60 I e y 



states with the degeneracy d, we express the exchange model 
(5.2) in terms of the d 2 - 1 independent operators X". With ex- 
tensive use of the generalized Stevens' multiplicative factors 
we visualize the wave functions using the formula (4.7). 
The static RPA susceptibility for the exchange model (5.14) 
is useful to determine the second-order phase transition line 
from the disorder phase under uniform external fields. The 
transition temperature is determined by (5.16). The GL free- 
energy expansion without the external fields (5.25) describes 
the entanglement of the multipoles in the ordered phase. The 
static RPA susceptibility in the ordered phase is given by 
(5.29). 

The physical multipole moments Q pq and M pq are evalu- 
ated by (3.17), in which the spherical tensor operators J pq are 
expressed as linear-combinations of the operators X" . The ex- 
istence of the multipole moments give rise to the electric and 
the magnetic fields near the magnetic ions, which are deter- 
mined by (2.12) and (2.13) with use of the vector spherical 
harmonics (A-l). The explicit example using Ce A Lai_ x B6 are 
given in §6. 

The analysis of the multipole exchange system tends to be 
complicated without systematic descriptions. The entangle- 
ment of the multipoles plays a key role to understand anoma- 
lous responses to external fields. The entanglement is com- 
posed concisely in the structure constants (5.4), which is use- 
ful to grasp a whole structure of the system. 

In recent years, the experimental techniques have been 
developed extensively to observe semi-quantitatively higher- 
rank multipoles. A quantitative analysis of a trace of the mul- 
tipole moments using the electromagnetic probes could accel- 
erate further development. 
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Appendix A: Vector Spherical Harmonics 

In contrast to a scalar field, the spatial rotation transforms 
not only the position r but also the direction of the vector 
field. The uniform vector is transformed as if it has the angular 
momentum ("spin") 1 . Consequently, it is natural to construct 
the vector spherical harmonics as a direct product of Y( m {r) 
and the spherical unit vector of rank 1, 



(6.25b) 



The observed thermodynamic anomalies in phase IV can be 
understood in this way. The detailed mean-field analysis 28 ' 
is consistent with the present discussions based on the GL 
expansion. 

7. Summary 

We have discussed the description of multipole degrees 
of freedom in the consecutive fashion. In the restricted CEF 



= z w 1 pq) y ^ (fi)ci 



= (-i)^ +1 V^TZ {- q I 

mm' ^ ' 

(A-l) 

where I = p, p ± 1 . It is easy to see by definition that Y ( pq (r) 
transforms like Y pq (r) under spatial rotation. Note that Y pq (r) 
is also an eigenfunction of the orbital angular momentum, 



t l Y c pq (,r) = ({(+\)Y e pq {r). 



(A-2) 
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If we introduce a "total" angular momentum from the orbital 
and the "spin" angular momenta, Y pq (f) is its eigenfunction, 
and the indices p and q represent the quantum numbers of the 
magnitude and the projection of the "total" angular momen- 
tum, respectively. 

The spherical unit vector is given by the cartesian unit vec- 
tors, 



/ Y p 

l, 1 pq 




e n = - 
ew = e 



-^=(e x + tey), 
V2 



Fig. A-l. The schematic relations between Y e pq , I and r. 



(A-3) 



e\-\ = —(e x - ie y ), 
V2 



-ifxY p+ q \f) = 



— £- — y p (f) 



(A-9c) 



which satisfy the orthogonality, e, , ■ e\ m - S mm >. This defi- „, , . . . . , „/ . , , 

3 ., , , 3 } m . . „ A The schematic relations between Y„„, t and r are shown 

nition is compatible with the real representation in §2.4. An ^ j , q 

arbitrary vector is expressed as 



Since the differential operator V also has a vector property, 
(A-4) derivatives can be expanded in terms of Y pq (f). With the help 
of the identity, V = f(d/dr) - (i/r)(f x (), we obtain 



where the spherical components of A are deduced from their 
cartesian components as similar to (A-3). The simplest cases 
of the vector spherical harmonics are 



y o ° (r) = 0, Yl{r) 



1 



1 



f. (A-5) 



1?v "~ VlT 19 ' 0,)V 

The complex conjugation is given by 

[Y%(r)J = {-\)™ +M Y ( p _ q {f). (A-6) 

The vector spherical harmonics are also expressed in terms 
of Y pq (f), f and t, 

1 



YUf) = 



1 pq 



(f) = 



Y p+ q \f) = 



1 

VP(2P + 1) 

-1 



V(p + m~P + 1) 



(Y pq (r), (A-7a) 
{pf-ifxt)Y pq (f\ (A-7b) 
[{p + l)f + if x <] F M (r). 



(A-7c) 



These expressions are derived from the definition, (A-l). 

To characterize the direction of YUf), the scalar and the 
vector products with f are useful. Using identities, f ■ I = 
f ■ (f x I) = and f x (r x () = -t, we obtain 



r ■ YUf) = 0, 



r-Y p ~\f) = 



r-Y p+ q \f) 



2p + l 



Y pq (f), 



P+l 



"\|2p + l 



YUf), 



(A-8a) 
(A-8b) 

(A-8c) 



and 



if X YUf) 



P + l Y p -\r) i 
^2p+\ pq {r) + 



2p + \ 



Y p :\f), 



-ifxY p -\f) = 



p+l 



<2 P+ i Yp ^> 



(A-9a) 
(A-9b) 







2p+ 1 





d 



+Y p+ q \f). 



'2p+ 1 \ r dr 



f(r) 
(A- 10) 



where f(r) is an arbitrary function of the radial coordinate r. 
The divergence and the rotation are also given by 



- V ■ (f(r)Y p pq (f)) = 0, 
-V-{f(r)Y p pq -\f)) = Y pq (f) 

-V-(f(r)Y p pq - l (f)) = Y pq (f). 
-iVx(f(r)Y p q (f)) = Y p q \f). 



(A- 11 a) 



Y p :\f) 



2p + l\ r 



N x I 





/p-1 


2/5+1 




p + l 1 


'p + 2 


2p + 1 ' 




p+l 




2/7+1 




>■ 


p+l 




V2p + 1 \ r 






2p + 1 





dr 



fir), 
(A- lib) 



+ d-rV^ 
(A- 11c) 



4m" 



(A- 12a) 



(A- 12b) 



f(r). 
(A- 12c) 



From those formula, we have some useful relations, 

V • (ffU) = V • {r p ~'Y p pq ) = V - = 0, (A13a) 

Vx(^^) = Vx[g] = 0, 



(A- 13b) 
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± pq 



p+l 



Vx 



y p 
1 pi 
irP +\ 



where j p {x) and h p \x) are the spherical Bessel function and 
the spherical Hankel function of the first kind, respectively, 
(A- 13c) and r< = min(r, r'), r > = max(r, r'). In the limit of k — > 0, the 
Green's function becomes 



V(p + D(2p+l) vJ>+1 

r p+2 r P1 ■ 



g{r,r') 



1 



(B-6) 



(A- 13d) 

A scalar product of two vector spherical harmonics are ex- 
panded in terms of Y( m (r) as 

Y e> * (r) ■ Y {l (r) = 



_ 1) „„ y (2pi-H)(2p 2 + l)(2f |t l)(2f;7 



1) 



2p + 1 ' 

Using the Green's function, we have the solution for (B-3) 

/p 
dr'^Z^pir'). (B-7) 
r '> 

For regions outside the source distribution, we have r< = r' 
and r> = r. Then, the multipole expansion for the scalar field 
is given by 







{ pi p 2 
■m -q\ q 2 



co p . 

^ = Z Z — z «wg«. ( B - g ) 

(A- 14) w ^ tn t ' le e l ectr i c multipole moment, 

Q pq = j drrPz; q (r)p{r). (B-9) 



Pi Pi 1 



Y tm (f\ 



where curly bracket is the 6j symbol. 13) Using this formula, 
we evaluate the angular integral, 



J drYUr)[Y%l(r)-Y% qi (r)\ = 



B.2 The vector potential 

The derivation here essentially follows that given by 
Schwartz in the context of the hyperfine structure of nuclear 



! j(2p! + V){2p2 + l)(2d + 1)(2( 2 + 1)(2( + 1) matter. 9 ' In the case of the gauge V • A = 0, the vector poten- 
i~ ' "Y 4^ tial is parallel to Y p pq (r) due to (A-ll). Namely, we express the 



I Pi 


Pi\ 


{pi 


h 


0- 


-m -qi 






Pi 


(A- 15) 



solution in the form, 



Mr) = X Z z «WfwW' 

P=0 q=-p 

In the case of I = m = 0, we have the projective orthogonality wne re we have introduced 
relation, 



dr Y e p \* qi (r) ■ Y l p \ qi (r) = S PlP2 6 qiqi 6 (l(l . 



(A-16) 



Zpq(f) = -i 



\p{2 P + \y 



tZ pq (r) 
ip 



(B-10) 



(B-ll) 



Appendix B: Details of Multipole Expansion 

B.l The scalar potential 

We write the solution of the Poisson equation in the form, 

oo p 

<P(r) = J]J]z pq (r)e pq (r). (B-1) 

p=0 q=-p 

Substituting this solution into the Poisson equation, we have 



for notational simplicity. Substituting this into the Poisson 
equation and using (A-2), we have 



pq 



y 2 PiP + 1) 



(pq(r) = -—j(r). (B-12) 



By the orthogonality, (A-16), we have the differential equation 

for £ pq (r), 



J]z pq (r) 



y2 P(P + 1) 



6 P qir) = -4np(r), (B-2) 



,2 PiP + 1) 



f«W = - 



p(2p + i) r 

TT7 J d ™pqir)-Kr). 



where is the radial part of the laplacian. Using the orthog- 
onality of the spherical harmonics, we obtain the differential 
equation for 9 pq (r), 



bpqv/ c(p+\) J "'~" ?v 

(B-13) 

We obtain the solution using the Green's function (B-6), 



/p 
dr'-^rZ*(r')-j(r'). (B 



14) 



y2 PiP + 1) 



pq (r) = -(2p + 1) 



\dfT pq 



(r)p(r). (B-3) 



In order to solve this equation, we consider the radial part of 
the inhomogeneous Helmholtz equation with a point source, 



7 2 + k 2_PiP^) 



(r,r') = - 



6(r - r') 



cip + 1) J 

For regions outside the source distribution, we have 

oo p . 

A(r) = Z Z —i Z p^ M p^ (BT5) 

P=0 9=-p 

where we have defined the magnetic multipole moment as 



The solution (the Green's function) is given by 
g {r,r') = ikj p (kr < )h p \kr > ), 



(B-4) 



(B-5) 



M 



pq 



cip 



°—JdrrP[z; q (f).j( r )]. (B 



16) 
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M pq = 4^Tp J drr"- l Y p pq u (r) ■ ^^)- (B-17) 

Now, we introduce the magnetization density by 

j(r) = cV x M(r). (B-18) 
With use of the identity, 

r -^- = rx(VxM) = V(r-M)-[l+(r-V)]M, (B-19) 



A(p,J,L,S) = (-1) J+L+S+P (2J + 1) 



\J J p 

[L L S 

~p L L 

L -L 



(C-3) 



Note that in the case of J = L + S , the coefficient becomes 
independent of L and S , 

-l 



A(p,J,L,S) 



and the integrations by part, we obtain 



p J J 
J -J 



(for 5 =J-L). (C-4) 



M(r) 



= J *v[r"z; ? (r)].M(r), (B-20) 

where we have used (A- 13a) and (A- 13c). 

The current density is expressed in terms of the angular- 
momentum density. The orbital current has the relation 



With the help of (C-2), the reduced matrix element of the 
electric multipole operator (3.2) is given by 

(nj\\Q p \\nj) = Y,(nJ\\Q P U)\\nJ) 

j 

= A(p, J, L,S)n (nLL | Qpoti) \ nLL) , (C-5) 



with 



QpoU) = -e drS(r - r 7 >% (r). 



(C-6) 



rxjorb.(r) = 2p B c{(r), 



while the spin current is / spin = 2yU B cV x s(r) by definition. 
The latter gives the spin magnetization density as M sp i n (r) = 
2p B s(r). Using (B-17) and (B-20) for the orbital and the spin 
parts, respectively, we obtain 



Here, we have used that (nLL | Qpo(j) \ nLL) is independent 
(B-21) of j, as will be shown shortly. 

Similarly, we apply the following relation to the magnetic 
multipole operator, (3.3), 



2t(r) 
p + l 



+ 2s(r) 



(B-22) 



M„=ita J drV[r p z; q (r)]- 

The orbital angular-momentum densities are expressed in 
terms of the /-electron operators as 

e( r ) = l^6(r- rj )e^ , S (r) = ^T6(r- rj ) S )i . 



(nJ J (f(r)Y { pq (r) ■ g(S)) p \\ nj} = (27 + 1) V^TI 

(J J p) 

x\L L (\(nL\\fY ( \\nL)(nS\\ gl \\nS) 



[S S 

= A(p,J,L,S,i) 



1 



tJ2p~TT 



iS(S + 1)(25 + 1) 
x(nLL\f(r)Y m (r)\nLL)(nS \\ gl \\nS) , (C-7) 



(B-23) where we have defined 
Appendix C: Derivation of Generalized Stevens' factors 

Let us express the wave function of the multiplet, J L2s+\, 

with the f configuration in the Russell-Sanders scheme as ^(p, J,L,S,i) = (27 + 1) y/S(S + 1)(2S + 1) 



\nJM) = (-!)— V57TT2(/m L m I) 



\nLm) \nScr) . 
(C-l) 



(J J p) 

J l l e\ 

If f l) 

t L L 

L -L) 
(C-8) 

13)1 



We first express the reduced matrix elements of the multipole and 3 x 3 curly bracket denotes the 9j symbol. li) Note that the 

operators in terms of the expectation value of the particular inner product Y e pq (r) ■ g transforms like Y pq (r) under spatial 

orbital state, 16 ' \nLL). Then, we derive the expectation value rotation. The result is 
for the Hund's-rule ground state with the maximum L. 



C.l The reduced matrix elements of multipole operators 

For an operator of rank p which acts only on the orbital part 
of the wave function, we have the relation, 6 ' 



(nj\\f p (L)\\nj) = 

(-1) J+L+S+P (2J + 1) j 
= A(p, J,L,S) (nLL | f p $ | nLLj , 



(nJ || M p || nJ) = Yj ( nJ J < b 0') + [M p p -\f) ■ Sj] p \\ 

j 

= A(p, J, L, S )n (nLL \ ( j) \ nLL) 

+ A(p,J,L,S,p-l)(nLL\f p - W (j)\nLL), (C-9) 



\J J 
\L L 



ji(nL||/„||nL) 



with 



(C-2) 



Mptif) = jfl J drd(r - r,)V (r"Z p0 (r)) ■ t h (C-10) 

f p - w (j) = 2p B ^p(4p 2 - 1) J" drd(r - rj )r»-%- W (r). 

(C-ll) 
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C.2 The expectation values of the Hund's-rule ground state 
The Hund's-rule ground multiplet in the Russell-Sanders 
scheme with /" configuration is characterized by the quantum 
numbers (JLS), 



K pk (n) = 7 



P 3 3\ 

o oj 



n-1 



7=1 v 



(k 3 3 
4-7 7-4 



J = \L-S\, L = J](4-7). 5 =n/2, 

7=1 



(C-12) 



(-iy 

- V7^ ] , (otherwise). (C-19b) 

From (C-5), (C-9) and (C-18), we have the reduced matrix 
elements of the multipole operators 



for n < 7, otherwise 



(nj\\Q p \\nj) 



J = L + S, L = Y^(4- j), S =7 -n/2. (C-13) 



7=1 



The orbital wave function is expressed by the Slater de- 
terminant of the one-body atomic wave function, y>4_ m (r) = 
Rf(r)Y 3m (r), 



-e (rP) 
(nj\\M p \\nj) {-\y +x p 
P 



(C-20) 



A(p, J,L,S)K pp (n), 
+ \ ^P(4p 2 ~ mp, J,L,S,p- l)^_ ljP _i(n). (C-21) 



\nLL) 



<^i(ri) <pi(r 2 ) ■■■ <fii(r„) 

<P2(ri) <P2(ri) ■■■ fi(r n ) 

<Pn(r\) f n ( r l) ■■■ <Pn(r„) 

<-\y>-i 



C.3 Generalized Stevens' factors 

Applying the Wigner-Eckart theorem to (3.1 1), we obtain 



y/n 



<Pm(rj)\n - l;mj), 



(C-14) 



where \n - l;mj) is (n - 1) X (n - 1) Slater determinant from 
which <p m and rj are eliminated. 

We consider an operator in the form, 

0(f) = J dr6(r - rj)f(r)g(r), (C-15) 

then we evaluate the expectation value as 
l(d(j))=n(nLL\d(f)\nLL) 

= </W>X f drY; m (r)g(r)Y 3m (r), (C-16) 

m ^ 

where (/(r)) is the radial average, and l(0( j)) is independent 
of j. Using the formula 



8 (p) = 

6n 



Sn 



(nj\\ Q p \\nj} 
-e(rP)(j\\j p \\jy 

(nj\\M p \\nJ~j 
p B (m-i) (J || Jp || 7)' 



(C-22) 



(C-23) 



With (C-20), (C-21) and (3.7), the generalized Stevens' multi- 
plicative factors gjf ) are evaluated in Table C-3 for the Hund's- 
rule ground multiplet with /" configuration, in which we also 
give the ratio of the orbital and the spin contributions to the 
magnetic multipoles, 

^(orbital) 



,00 



^(spin) 



(C-24) 



Note that the orbital contribution in rj, 3) vanishes for L = 5 
due to Kxiiri) = 0. 

In the case of p = 1, we have 

L 

K QQ (n) = n, K Ql (n) 



! 



dfY pigi (f)Y P2q2 (f)Y P3q3 (f) = 



(2 Pl + 1)(2 P2 + 1)(2 P3 + 1) 



An 



Co' 


Pi 


Pl\ 




Pi 


pi\ 











qi 


?3/ 



(C-17) 



and (A-15), we obtain the required expectation values, 
l{Qpo(jj) = -e{r p )Kpp(n), (C-18a) 



A(1,J,L,S) _ 1 
^(1,7,L,5,0) 



1 + 



2V3' 
L(L+ l)-5(5 + 1) 



1 



and we obtain 



1 



7(7 + 1) 
L(L+ l)-S(S + 1) 



2V3[ ■/(•/+!) 
3 L(L+ l)-S(S + 1) 



(C-25) 



(C-26) 



= /"b (^')(-ir'^/vi»- (c-isb) 



/(/t>-io(7)) = A<b (r p - 1 )2 y fp~(4p^)Kp- l p- l (n), (C-18c) 
where we have defined 

w= 7 (s o)g(o /-j ji 4 y-v. 



(C-19a) 



*m = -- 

Sn 2 7(7 + 1) 

This is nothing but the Lande's g factor. 
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Table C I. The tesseral harmonics (multiplied by r p ) in the cartesian coordinate. The parity is given by (—l) p . 



p q rPZf q {r) or r"Z p0 (r) 



1 



1 z 
1 x 



3 -z(5z 2 -3r 2 ) 
2 ^-/) 



2 

VTo 



3v 



4 | (35 Z 4 - 30z 2 r 2 + 3r 4 ) 

2 x(^ 2 )( 7 ^ 2 ) 

" ^zxfx 2 ^ 2 ) 



4 
V35 



(/- 6 xV +/) 



5 - (63z 5 - 70z 3 r 2 + 15zr 4 ) 

O J 

^px[r 4 + 7z 2 (3z 2 -2r 2 )] 
z(x 2 -y 2 )(3z 2 -r 2 ) 



1 

2 
3 
4 
5 



O05 



^(, 2 -3/)(9z 2 -r 2 ) 

^z(x 4 -6xV + /) 
3Vl4 



06- 



x(x 4 -10xV + 5/) 



2 ^(3z 2 -r 2 ) 




1 V3zx 


V3yz 


2 f 


Vlry 



4 

VHryz 

VTo 



-yz 



(7z 2 -3r 2 ) 



VTo 

f^-r 2 ) 

Vvo , 2 2 x 
— Jz(3x -y) 

V35 



PS) 



VT5 



005 



6 — (231z 6 - 315z 4 r 2 + 105z 2 r 4 - 5r 6 ) 
16 v 1 

^zx[5r 4 + 3z 2 (llz 2 -10r 2 )] 



2 
V70 

3_V35 
2 

3 VTi 
UK. 



y[r 4 + 7z 2 (3z 2 -2r 2 )] 
'■xyz(3z 2 - r 2 ) 
(3x 2 -f)(9z 2 -S) 
xyz(x 2 -y 2 ) 
-y(5x 4 -\0x 2 y 2 + /) 



1 

2 
3 
4 
5 
6 



V2T 



yz[5r 4 + 3z 2 (llz 2 -10r 2 )] 



V210 
32 

^zx(x 2 -3/)(llz 2 -3r 2 ) 

1^1^(^-10^ + 5/) 
16 v ' 

^[/-15xV(x 2 -/)-/] 



( x 2 _ ^ [ r 4 + 3z 2 ^ U 2 _ 6r 2j] .V210 xy j r 4 + ^ ^2 _ ^ 

yz(3x 2 -/)(llz 2 -3r 2 ) 
xy(x 2 -y 2 )(Uz 2 -r 2 ) 
yz(5x 4 -10xV+/) 
xy(3x 4 - lOx 2 / + 3/) 
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Table C-2. The cubic harmonics as linear combinations of the tesseral harmonics. The parity is given by (-l) p . 
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Table C-3. The generalized Stevens' factors for the Hund's-rule ground multiplet in the Russell-Saunders scheme. 
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